home *** CD-ROM | disk | FTP | other *** search
- /**********************************************************************
- Algorithm for getting 3-D coordinates
- from 2-D coordinates in mirror symmetric model
-
- Copyright (C) 1993 Shusaku Furushima.
-
- Caution: This source code cannot be executed itself. Please understand
- the concepts of this algorithm and modify this source code to fit your
- products.
-
- This source code is copyrighted, but you can freely use all or a part of
- source code for your products (commercial use also allowed).
- **********************************************************************/
-
- #include <stdio.h>
- #include <math.h>
-
- #define real double
- /* max point of qr analysis */
- #define RQRMP 10
- #define RQRNP 30
-
- /******** Rcoord ********/
- typedef struct {
- real x;
- real y;
- real z;
- } Rcoord;
-
- /******** 4 x 4 matrix **********/
- typedef real Rmat[4][4];
-
- /********** mathematical functions ********/
- extern void Rqranalysis(int m, int n, real qa[][RQRNP],
- real qy[], real qz[]);
-
- /******** general 4x4 matrix *********/
- class Rmat44{
- void initmat(Rmat amat);
- public:
- void matmul(Rmat mat2);
- Rmat mat;
- void IRmat44();
- void rotp(real p);
- void rotq(real p);
- void rotr(real p);
- void copyfrom(Rmat mat2);
- void move(real x, real y, real z);
- void bepers();
- void mul(real x, real y, real z);
- void getcoord2(Rcoord *coord, Rcoord *cd, int pers);
- };
-
- /********* 4x4 matrix for display objects *********/
- class Rmatobj{
- public:
- Rmat44 *devmat; /* matrix of [3D -> 2D screen] */
- Rmat44 *symm; /* symmetrical matrix of devmat */
- real pp, pq, pr; /* rotate angle */
- real dx, dy, dz; /* move coordinate */
- real bx, by; /* scale to screen 2d coord */
- int pers; /* flag if pers:1 or ortho:0 */
-
- void IRmatobj(); /* init matrixes */
- void makedev(); /* make devmat */
- void setsymmatx(); /* set symm matrix of devmat */
- Rcoord getdev(Rcoord *pt); /* get screen coord */
- Rcoord calc3d(Rcoord *pt1, Rcoord *pt2); /* get 3D coord from 2D coord */
- };
-
-
-
- /********** least square method by QR (without pivoting) ***********/
- /********* qranalysis simple version ********/
- /* m number of unknown param. */
- /* n number of data sets. */
- /* qa data param matrix. */
- /* qy data value vector. */
- /* qz returning answer. */
- /* qrnp define value of points number. */
- /* qrmp define value of param number. */
- /******************************************/
- void Rqranalysis(int m,int n,real qa[][RQRNP],
- real qy[],real qz[])
- {
- int i,j,k;
- static real a,qq[RQRMP][RQRNP],rq[RQRMP][RQRMP];
-
- for( k = 1; k <= m; k++)
- {
- a = 0;
- for(i = 1; i <=n; i++)
- a = a + qa[k][i] * qa[k][i];
- rq[k][k] = sqrt(a);
- for( i = 1; i <= n; i++)
- qq[k][i] = qa[k][i] / rq[k][k];
- for( j = k+1; j <= m; j++)
- {
- a = 0;
- for( i = 1; i <= n; i++)
- a = a + qq[k][i] * qa[j][i];
- rq[j][k] = a;
- for( i = 1; i <= n; i++)
- qa[j][i] = qa[j][i] - qq[k][i] * rq[j][k];
- }
- }
-
- for( i = 1; i <= m; i++)
- {
- a = 0;
- for( j = 1; j <= n; j++)
- a = a + qq[i][j] * qy[j];
- qz[i] = a;
- }
- qz[m] = qz[m] / rq[m][m];
- for( i = m - 1; i >= 1; i--)
- {
- for( j = i + 1; j <= m; j++)
- qz[i] = qz[i] - qz[j] * rq[j][i];
- qz[i] = qz[i] / rq[i][i];
- }
- }
-
- /******** mathod for general matrix Rmat44 *********/
- /*** initialise matrix ***/
- void Rmat44::IRmat44()
- {
- initmat(mat);
- }
-
- /*** set [I] matrix ***/
- void Rmat44::initmat(Rmat amat)
- {
- int i, j;
- for (i = 0; i < 4; i++){
- for (j = 0; j < 4; j++){
- amat[i][j] = 0.0;
- if (i == j)
- amat[i][j] = 1.0;
- }
- }
- }
-
- /*** rotate around Y axis ***/
- void Rmat44::rotp(real p)
- {
- Rmat tmp;
-
- initmat(tmp);
- tmp[0][0] = cos(p);
- tmp[2][0] = sin(p);
- tmp[0][2] = - sin(p);
- tmp[2][2] = cos(p);
- matmul(tmp);
- }
-
- /*** rotate around X axis ***/
- void Rmat44::rotq(real p)
- {
- Rmat tmp;
-
- initmat(tmp);
- tmp[1][1] = cos(p);
- tmp[2][1] = -sin(p);
- tmp[1][2] = sin(p);
- tmp[2][2] = cos(p);
- matmul(tmp);
- }
-
- /*** rotate around Z axis ***/
- void Rmat44::rotr(real p)
- {
- Rmat tmp;
-
- initmat(tmp);
- tmp[0][0] = cos(p);
- tmp[1][0] = -sin(p);
- tmp[0][1] = sin(p);
- tmp[1][1] = cos(p);
- matmul(tmp);
- }
-
- /*** shift coord ***/
- void Rmat44::move(real x, real y, real z)
- {
- Rmat tmp;
-
- initmat(tmp);
- tmp[0][3] = x;
- tmp[1][3] = y;
- tmp[2][3] = z;
- matmul(tmp);
- }
-
- /*** set to perspective matrix ***/
- void Rmat44::bepers()
- {
- Rmat tmp;
-
- initmat(tmp);
- /* project plane is always z=1 and eye point is x=y=z=0 */
- tmp[3][2] = 1.0;
- matmul(tmp);
- }
-
- /*** magnify value ***/
- void Rmat44::mul(real x, real y, real z)
- {
- Rmat tmp;
-
- initmat(tmp);
- tmp[0][0] = x;
- tmp[1][1] = y;
- tmp[2][2] = z;
- matmul(tmp);
- }
-
- /*** multiple two matrices ***/
- void Rmat44::matmul(Rmat mat2)
- {
- int i, j, k;
- Rmat tmp;
-
- for (i = 0; i < 4; i++){
- for (j = 0; j < 4; j++){
- tmp[i][j] = mat[i][j];
- }
- }
- for (i = 0; i < 4; i++){
- for (j = 0; j < 4; j++){
- mat[i][j] = 0.0;
- for (k = 0; k < 4; k++){
- mat[i][j] += tmp[k][j] * mat2[i][k];
- }
- }
- }
- }
-
- /*** copy matrix ***/
- void Rmat44::copyfrom(Rmat mat2)
- {
- int i,j;
- for (i = 0; i < 4; i++){
- for (j = 0; j < 4; j++){
- mat[i][j] = mat2[i][j];
- }
- }
- }
-
- /*** get 2D coordinates from 3D coordinates from the matrix ***/
- void Rmat44::getcoord2(Rcoord *coord, Rcoord *cd, int pers)
- {
- int i,j;
- real w;
-
- coord->x = mat[0][0] * cd->x + mat[0][1] * cd->y +
- mat[0][2] * cd->z + mat[0][3];
- coord->y = mat[1][0] * cd->x + mat[1][1] * cd->y +
- mat[1][2] * cd->z + mat[1][3];
- coord->z = 0.0;
- if (pers){
- w = mat[3][0] * cd->x + mat[3][1] * cd->y +
- mat[3][2] * cd->z + mat[3][3];
- w = 1.0 / w;
- coord->x *= w;
- coord->y *= w;
- }
- }
-
- /************* method for Roku matrix Rmatobj **************/
- /*** initialize matrices ***/
- void Rmatobj::IRmatobj()
- {
- devmat = new(Rmat44);
- symm = new(Rmat44);
- devmat->IRmat44();
- symm->IRmat44();
- /*** You can change these values in your application ***/
- pp = -0.9; /* default angle */
- pq = -0.25;
- pr = 0.0;
-
- dx = 0.0; /* default position */
- dy = 0.0;
- dz = 1.5;
-
- bx = by = 100.0; /* default magnify value */
-
- pers = 1; /* default is perspective (not orthogonal) */
- }
-
- /*** make matrix of [3D coordinates -> 2D screen coordinates] ***/
- void Rmatobj::makedev()
- {
- devmat->IRmat44();
- devmat->rotp(pp); /* rotate around Y axis */
- devmat->rotq(pq); /* rotate around X axis */
- devmat->rotr(pr); /* rotate around Z axis */
- devmat->move(dx, dy, dz); /* shift position */
- if (pers)
- devmat->bepers(); /* set to perspective matrix */
- dsmat->mul(bx, by, 1.0); /* magnify to be screen scale */
- }
-
- /*** get 2D screen coordinates from 3D coordinates ***/
- Rcoord Rmatobj::getdev(Rcoord *pt)
- {
- Rcoord tmp;
-
- devmat->getcoord2 (&tmp, pt, pers);
- return(tmp);
- }
-
- /*** set plane symmetry matrix ***/
- void Rmatobj::setsymmatx()
- {
- symm->IRmat44();
- symm->mat[0][0] = - 1.0;
- symm->matmul(devmat->mat);
- }
-
- /*** calc 3D coordinates from two 2D coordinates ***/
- Rcoord Rmatobj::calc3d(Rcoord *pt1, Rcoord *pt2)
- {
- Rcoord p;
- real qa[RQRMP][RQRNP], qy[RQRNP], qz[RQRMP];
-
- qa[1][1] = devmat->mat[0][0] - devmat->mat[3][0] * pt1->x;
- qa[2][1] = devmat->mat[0][1] - devmat->mat[3][1] * pt1->x;
- qa[3][1] = devmat->mat[0][2] - devmat->mat[3][2] * pt1->x;
- qy[1] = - devmat->mat[0][3] + devmat->mat[3][3] * pt1->x;
-
- qa[1][2] = devmat->mat[1][0] - devmat->mat[3][0] * pt1->y;
- qa[2][2] = devmat->mat[1][1] - devmat->mat[3][1] * pt1->y;
- qa[3][2] = devmat->mat[1][2] - devmat->mat[3][2] * pt1->y;
- qy[2] = - devmat->mat[1][3] + devmat->mat[3][3] * pt1->y;
-
- qa[1][3] = symm->mat[0][0] - symm->mat[3][0] * pt2->x;
- qa[2][3] = symm->mat[0][1] - symm->mat[3][1] * pt2->x;
- qa[3][3] = symm->mat[0][2] - symm->mat[3][2] * pt2->x;
- qy[3] = - symm->mat[0][3] + symm->mat[3][3] * pt2->x;
-
- qa[1][4] = symm->mat[1][0] - symm->mat[3][0] * pt2->y;
- qa[2][4] = symm->mat[1][1] - symm->mat[3][1] * pt2->y;
- qa[3][4] = symm->mat[1][2] - symm->mat[3][2] * pt2->y;
- qy[4] = - symm->mat[1][3] + symm->mat[3][3] * pt2->y;
-
- Rqranalysis(3, 4, qa, qy, qz);
- p.x = qz[1];
- p.y = qz[2];
- p.z = qz[3];
- return(p);
- }
-
- /************** summary of 2D-3D conversion routine in HexaGreat-3D *****************/
- void main()
- {
- Rmatobj *rmat;
-
- Rcoord gridscpt[MAX], gridpt[MAX]; /* Grid data */
- int ngpoint, ngline, gridline[MAX][2];
-
- Rcoord inputscpt[MAX][2], inputpt[MAX][2]; /* Drawing data */
- int ninputpt, ninputln, inputline[MAX][2];
-
- short x, y;
- int i;
-
- :
- // set grid coordinates and line informations.
- :
- rmat = new(Rmatobj);
- rmat->IRmatobj(); /* initialize matrix */
- :
- // changing display parameters (such as rmat->pp)
- :
- rmat->makedev(); /* get conversion matrix */
-
- for (i = 0; i < ngpoint; i++){ /* calculate screen coordinates of grid */
- gridscpt[i] = rmat->getdev(&gridpt[i]);
- }
- for (i = 0; i < ngline; i++){ /* display grid to screen */
- x = (short)(gridscpt[ gridline[i][0] ].x);
- y = (short)(gridscpt[ gridline[i][0] ].y);
- MoveTo(x, y);
- x = (short)(gridscpt[ gridline[i][1] ].x);
- y = (short)(gridscpt[ gridline[i][1] ].y);
- LineTo(x, y);
- }
- :
- // input 2D drawings.
- // at this moment, inputscpt[i][0] contains 2D screen coordinate of
- // left side point, and inputscpt[i][1] contains 2D screen coordinate
- // of right side point.
- // (only x and y values are used. z value is not used.
- :
- rmat->setsymmatx(); /* set symmetric matrix */
- for (i = 0; i < ninputpt; i++){
- /* calculate 3D coordinate */
- inputpt[i][0] = rmat->calc3d(&inputscpt[i][0], &inputscpt[i][1]);
- /* invert X axis sign */
- inputpt[i][1] = inputpt[i][0];
- inputpt[i][1].x = - inputpt[i][1].x;
- }
- :
- // consequently, inputpt[i][0] contains 3D coordinate of left side point,
- // and inputpt[i][1] contains 3D coordinate of right side point.
-
- }
-
-
-
-